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ABSTRACT 

We have performed a WISE (Wide-Field Infrared Survey Explorer) based study to iden¬ 
tify and characterize young stellar objects (YSOs) in 12° x 12° Perseus OB 2 association. 
Spectral energy distribution (SED) slope in range of 3.4 — 12/tm and a 5cr selection criteria 
were used to select our initial sample. Further manual inspection reduced our final catalog 
to 156 known and 119 YSO candidate. The spatial distribution of newly found YSOs all 
over the field shows an older generation of star formation which most of its massive mem¬ 
bers have evolved into main sequence stars. In contrast, the majority of younger members 
lie within the Perseus molecular cloud and currently active star forming clusters such as 
NGC 1333 and IC 348. We also identified additional 66 point sources which passed YSO 
selection criteria but are likely AGB stars. However their spatial distribution suggests 
that they may contain a fraction of the YSOs. Comparing our results with the commonly 
used color-color selections, we found that while color selection method fails in picking up 
bright but evolved weak disks, our SED fitting method can identify such sources, including 
transitional disks. In addition we have less contamination with background sources such 
as galaxies, but in a price of loosing fainter ( Jmag > 12) YSOs. Finally we employed a 
Bayesian Monte Carlo SED fitting method to determine the characteristics of each YSO 
candidate. Distribution of SED slopes and model driven age and mass confirms separated 
YSO populations with suggested three age groups of younger than 1 Myr old, 1 — 5 Myr 
old, and older than 5 Myrs which agrees with the age of Per OB 2 association and currently 
star forming sites within the cloud. 


INTRODUCTION 


Perseus is one of the most interesting and well studied star forming complexes, containing massive 
stars (Per OB2 association), young stellar clusters (e.g. IC348 and NGC1333) and several dense star 
forming cores within Perseus molecular cloud (e.g. |Kirk et ah ( 2006| ) ; [Ridge et ah] (2006); Rosolowsky 
et al. (2008)). Perseus OB2 is the second closest OB association to the Sun (after Scorpius-Centaurus 


association) at a distance of ^ 300 pc. de Zeeuw et al. (1999) estimated that Per OB2 is approximately 
6 Myr old while other studies also suggest different generation of star formation with an age less than 


15 Myr (e.g. Bally et al. (2008) and references therein). 


Per OB2 contains both massive and intermediate mass stars while low mass to intermediate mass 


stars are currently forming within the Perseus molecular cloud (e.g. Enoch et al. (2009); Arce et al. 


(2010); Sadavoy et al. (2014])). Therefore identifying the members of the association and determining 


1 Harvard-Smithsonian Center for Astrophysics, Cambridge, MA 02138 











































- 2 - 


their characteristics such as their age and mass will provide critical information on the process of 
different generations of star formation in this region. 

We constructed a panchromatic census of young stellar objects (YSOs) in Per OB2 using Wide- 
field Infrared Survey Explorer ( WISE) point source catalogue. We also employed all available optical 
and infrared all sky surveys including Spitzer , 2MASS and SDSS to construct spectral energy dis¬ 
tribution (SED) of each source in order to determine its physical parameters. We have performed a 
WISE based survey to identify young stellar objects throughout a 12° x 12° area covering the Per 
OB2 association. 


The Perseus molecular cloud covers a 2° x 6° 10 x 30 pc) area and was observed with the Spitzer 

space observatory in both IRAC (3.6 — 8.0 /un) and MIPS (24 — 160 pm) as part of “From Molecular 
Cores to Planet-Forming Disks” (c2d) Legacy Survey ( Jprgensen et al.|[2006 ; Rebuff et al. 2007). We 
selected the c2d covered region as footprint, but did not limit our survey to the c2d area. c2d only 
covers very young star forming regions such as IC348 and NGC1333 clusters. To study the complete 
population of Perseus OB2 association we need to track all the members of the original ensemble, in 
particular the entire OB association which has blown up the initial molecular cloud and other young 
stellar objects, now spread into a larger area. Assuming an average speed of 1 — 2 km s -1 for stars, 
members can travel more than 20 pc in 10 Myrs in each direction which converts to an extended area 
of ^ 10° from the cloud center at 250 — 300 pc distance. 


To be complete, we selected a 12° x 12° region (50° <RA< 64° and 26° <Dec< 38°) around 
the Perseus c2d footprint which completely covers the Perseus OB2 association. This region slightly 
overlaps with Taurus star forming region in South-West and may contain older members of Pleiades 
in South as well. Figure [l] shows the c2d footprint, and the area selected around it in our survey. 
The Blue dots present WISE point sources. Unfortunately we do not have the spectral or dynamical 
data for each source, hence cannot confirm if each individual source in our final catalog is associated 
with the Per OB2 or the molecular cloud, but the luminosity and other parameters matches with 
the properties we expect for YSOs at the estimated distances. In addition, we are looking out of 
the Galactic plane and no other known background star forming region is along the line of sight. 
Foreground young stars, if any, would be much brighter than our criteria. 

Section [2] explains the data collection and data reduction procedures. In section [4] we describe 
our YSO candidate selection method and compare the results with other YSO identification schemes. 
Section [4] contains the analysis of our sample of YSO candidates. In Section [5] we describe our SED 
fitting modeling and the results. We summarize this study in Section |6j 


2. Archival Data and Data Reduction 

To identify young stars we established a multi-wavelength database of point sources in Per OB2 
star forming region, using all available large sky surveys. New YSO candidates then were selected 
through evidence for infrared excess in their spectral energy distributions (SED); specifically, we used 
power-law fits to the longer infrared wavelengths covered by the WISE all sky survey. In section [4] we 
detail the candidate selection procedure. 
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To summarize, WISE is an all sky survey in four infrared wavebands 3.4, 4.6, 12, and 22 fim 
(hereafter Wl, W2, W3 and W4) with an angular resolution of 6-6.5" for W1 to W3 and 12. 0" for W4. 
The WISE depth is not homogenous, but achieved 5cr point source sensitivities of 0.068, 0.098, 0.86 
and 5.4 mJy (16.6, 15.6, 11.3, 8.0 Vega mag) in unconfused regions in the ecliptic plane for the four 
bands respectively. Our point source database was initially constructed using the WISE preliminary 


data release (Wright et al. 2010), but we upgraded to the final source release in March 2012. The 
12° x 12° selected region in our census of Per OB2 includes more than 1.6 million WISE point sources 


(Figure [I]). 


After an initial inspection we realized that W4 sensitivity is very low for faint sources and in 
particular any W4 profile magnitude >8 should be considered with caution. We therefore limited our 
search for infrared excess to the range of 3 — 12 /mi, and restricted our source list to have S/N>7 in 
the three corresponding Wl, W2, and W3 bands. This criteria decreased the number of candidates to 
65,662. We did not use the WISE catalog flags to further filter our candidate list because we found 
the flags to be unreliable; in particular many sources were not flagged as extended or questionable 
when in fact we found that no point source actually existed. Instead all of the final candidates were 
visually inspected to remove contaminated or false detections. 


Our primary goal was the building of composite spectral energy distributions (SED) using archival 
data products. We initially considered WISE as the base catalog to astrometrically match with other 
data sets, but the WISE beam is significantly larger than other relevant all sky surveys. We also 
decided to use a near-IR band as a tracer of the stellar spectral type. Therefore we adopted the 
2MASS (Skrutskie et al. 2006) catalog and coordinates as the base catalog. Requiring our sources to 
have a 2MASS companion reduced the size of our sample to 48,692. 


2.1. Additional ancillary data 


The Perseus star forming region is well covered by many large sky surveys (see Figure |T]) from 
infrared to X-ray including: 2MASS, UKIRT Infrared Deep Sky Survey (UKIDSS), Spitzer (e.g., c2d, 
Evans et al. ( |2003D ), the InfraRed Imaging Surveyor (AKARI) ( |Murakami et al.||2007| ), Sloan Digital 


Sky Survey (SDSS, DR8+), the AAVSO Photometric All-Sky Survey (APASS, Release 7)Q Chandra 
X-Ray Observatory and X-ray Multi-Mirror Mission (XMM-Newton) as well as catalogs based on older 
scanned photographic plates (US Naval Observatory: USNO-B1, UCAC3) and the USNO-B1/2MASS 
re-reduction, PPMXL (Roeser et al. 2010). 


SDSS and APASS catalog measurements are adopted at optical wavelengths if available; otherwise, 
we have used PPMXL cataloged photometry for B and R bands. Figure [2] presents APASS and SDSS 
coverage of our surveyed region, showing how APASS provided nearly complete, deep optical coverage. 
While PPMXL covers the entire area, it is intrinsically based on USNO-bl and contains much larger 
errors than APASS. Because PPMXL contains corrected positions, proper motions and optical-infrared 
matches of stars, we used it to break degeneracies for multiple matches between different catalogs. 
X-ray observations were very concentrated toward the known star forming regions, and did not cover 
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the entire region. We benefited from archival X-ray data in publications to identify known YSOs. 
Because of its limited spatial scope (Figure [l]) we neglected the UKIDSS catalog entirely. AKARI flux 
measurements were consistent with IRAC and MIPS observations, but only a few of our candidates 
had high quality AKARI data, therefore this catalog was not used in final SED fitting. 


There are several narrow field catalogs at different wavebands in Perseus region, but to keep a high 
degree of uniformity in source detection and flux calibration we do not account for such observations. 
We further downloaded all SIMBAD sources in the survey area to identify previously known YSOs or 
census contaminating sources such as planetary nebulae, background galaxies or other type of known 
evolved stars. 


We employed the CDS X-match Service ] to match different catalogs and investigate the astro- 
metric uncertainties. This service employs the coordinate errors to calculate the probability density 
of two sources from two different catalogs. The probability density is given by convolution between 
two Gaussian distributions around each source Pineau et al. (201l]). We consider a completeness of 
99.7 % at 3 <t criterion (only 0.3% of the counterparts could be missed). 


3. IDENTIFYING DISK CANDIDATES 


We adopt the shape of the spectral energy distribution (SED) to identify sources with excess 
infrared emission as these are probable young stars with passively reradiating circumstellar disks 


(Lada 1987). Following on the Spitzer based SED analysis of Lada et al. (2006) and Muench et al. 


(2007), we used a simple power-law, least-squares fit of the 3 — 12 /am portion of the YSO’s SEDs as 


traced by the Wl, W2, and W3 WISE bands (3.4, 4.6 and 12 /im, respectively). 


While Muench et al. (2007) showed that the Spitzer 3 — 8 fi m index was relatively insensitive 


to extinction, we acknowledge that including 12 microns into the SED slope fit will result in larger 
intrinsic dispersions in the slope distribution because of variations in the silicate feature. Unreddened 
classical T-Tauri stars with SiO 4 in emission will have very large slopes. We do not anticipate silicate 
absorption due to Ay to “remove” a source from our sample, by lowering the measured 3 — 12 micron 
slope. 

Figure [3] presents a comparison between IRAC ass and WISE 0 ^ 3-12 for 135 sources which have 
both data set in our sample. Solid line shows a least square fit (with slope of 0.92 ± 0 . 2 , excluding 
three data points at right bottom corner with poor IRAC photometry) and dashed line presents equal 
values. There is a large scatter due to different wavelength, band width for each wavelength, beam size 
and sensitivity for Spitzer and WISE. A selection bias is involved as well. In Figure [3] all IRAC data is 
included if available, but only WISE with SNR>7 sources are selected in our initial sample. In general 
for ass < —2, we have found larger slopes for WISE 0 ^ 3 - 12 , where anemic disks and stars are located 


(discussed further in section 3.3). IRAC and WISE slopes are not expected to match completely, for 
this reason we are not reporting classification for individual objects based on categories defined in 


Muench et al. (2007). 


2 http://cdsxmatch.u-strasbg.fr/xmatch, Strasbourg Astronomical Data Center 
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The 22 /im was not considered in the initial slope-fit for two reasons: a) to keep detection 
uniformity for sources which have not been detected or have very low SNR in W4 and b) to avoid 
the overestimation of slopes due to poor sensitivity of W4 for faint sources with W4 > 8. However 
22 /im flux, when regarded as real, is considered in final complete SED modeling to study the disk 
characteristics. 

In the first attempt we filtered candidates based on WISE flags. We found a few strong YSO 
candidates which were located nearby saturated stars or within nebulous regions, had been flagged as 
poor detection or photometry due to contamination. On the other hand in visual inspection we found 
sources with very poor detection in cloudy areas which had high S/N fluxes reported with WISE high 
photometry quality flags but they did not survive our visual or point source profile inspections. Finally 
we decided to ignore the WISE flags at this stage and reconsider them after manual inspection. 

We considered all sources with S/N>7 in first 3 WISE bands in our candidate pool and fitted a 
power-law to each source. Figure [4] shows the calculated slopes versus 2MASS J magnitude for 48,692 
sources with S/N>7 in first 3 WISE bands in grey dots. To select the sources with excess emission 
we looked at the slope distributions in bins of 0.5 J magnitude and fitted a Gaussian to the slopes 
distribution within each interval. J band is the shortest wavelength that we have uniform data for, 
and is thus the least contaminated wavelength by circumstellar disk emission, which makes it the best 
tracer of photospheric luminosity. 

The selection criteria can be considered to be a convolution of variation in the intrinsic (photo- 
spheric) color as a function of apparent magnitude and the WISE sensitivity function. We did not 
consider disentangling these two functions as a useful exercise, and instead posited all sources with an 
excess greater than 5 sigma from the typical color as evidence for infrared excess. 

The black locus in Figure [4] shows a 5cr limit above the peak of distribution for each bin. All 
sources above the locus have significant excess in slope (and hence in IR emission) and are considered 
as candidates with excess. On average the slope increases with the Jmag after magnitude 12. Fainter 
objects in J band are either embedded sources, or background galaxies. Hence, faint YSOs are missed 
in this selection, but in return the final sample is less contaminated with background galaxies. 

The Jmag error is not considered in calculating the locus, but it will not dramatically change the 
results. The maximum error in Jmag in our entire data set is only 1.7%. The locus has been drawn 
based on 5a above the mean population in each J m ag bin which means a confidence of 99.999943% 
of IR excess. Considering the largest error for J m ag around the locus, selected candidates are still 
more than 4a above the mean, which translates into a 99.994% confidence in IR excess. Therefore 
by neglecting J m ag error in our selection procedure we are not maintaining the 5a confidence for all 
selected YSO candidates, but all of the selected candidates have statistically significant IR excess. 

669 sources met this criteria. However we are aware that the sample may contain none-stellar 
objects, contamination and false detections. In next step we filter candidate by manually inspecting 
them in optical and infrared images. 
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3.1. Visual inspection 


WISE and 2MASS and if available SDSS images of each candidate which passed the locus selection 
were manually inspected. Visual inspection revealed that half of the selected sources cannot be 
considered resolved point sources. The first group were rejected as being identified as extended objects 
like galaxies and the second portion could not be identified as distinguished point sources in nebulous 
regions in a combination of at least any three WISE bands. We performed the eye inspection in 
three rounds and categorized the sources as contaminated, saturated, extended and not resolved. 
All extended objects were identified as galaxies and PNe in further investigations and were rejected. 
Candidates that could not be resolved as point sources in at least any three WISE bands needed to be 
inspected for the luminosity profile. Contaminated, saturated and not resolved sources were accepted 
and back to the list if they were resolved in 2MASS and if luminosity profile presented a point source 
in at least three WISE bands. 


354 out of 669 candidates survived the visual inspection. These final candidates are presented 
by open red circles (new candidates) and blue triangles (known candidates) in Figure |4j The known 
AGB stars are presented with black dots. The sources in our list which are less likely to be YSOs 
and closer to the AGB stars are noted by black crosses. They are separated from the YSO candidates 
employing magnitude diagrams and will be discussed in more details in section |4j The known OB 


stars are also shown as yellow dots. The only OB star which is above the locus is BD + 30540. Rebull 


et al. (2007) identified this YSO as a class I source in their Spitzer c2d survey. It is a B8V star and 


have an estimated mass of 2.6 Mq. 


Figure [5] presents the spatial distribution of 354 final candidates on a WISE multi-color image. 
The brightness of each source is proportional to the the Wl-3 slope; sources with larger slopes, which 
presumably have stronger disks, appear brighter in Figure [5j 


3.2. Comparing classification schemes 


Rebull et al. (2011) identified new YSO candidates in Taurus-Auriga region using preliminary 


WISE data. Nearly 25% of our survey area is covered in their study, therefore we looked for common 


candidates to evaluate our detection method accuracy. Rebull et al. (2011) report 94 new YSO can¬ 


didates in Taurus-Auriga region, 22 of which lie within our survey area. Astrometry and photometry 
of more than 25% of sources have changed in WISE all-sky release, especially for bright sources in 
our list. Therefore preliminary and all sky release WISE catalogs do not contain identical IDs for 
all sources. Instead of matching IDs we performed a position match. Some photometries in the new 
WISE catalog also have sightly changed, but they are consistent enough with the preliminarily data to 
help selecting the right match. Rebull et al. catalog is plotted over our selection diagram in Figure [6j 
Our complete catalog of 354 sources are shown in black triangles. Rebull et al. known YSOs are shown 
with green dots and their new YSO candidates are shown with blue dots. They have also rejected 
many sources (shown in magenta dots) after various inspections. We have identified only 8 of Rebull’s 
new YSO candidates in our field. 14 unidentified sources are almost fainter than 12 magnitude in 
2MASS J band and despite color excess they lie below our selection criteria. 
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We also have 17 of Rebull’s known candidates in our field. 14 have been classified as YSO 
candidates in our sample, one does not have 2MASS match and therefore not selected in our sample 
and two are below our selection criteria. 

Rebull et al. also have 129 rejected candidates in our field which passed their color cuts in first 
place. 82 of their rejected sources (as background galaxies) match with our initial sample but lie below 
our selection criteria. 35 do not have 2MASS point source match and therefore not selected in our list. 
We have 8 common sources in our visually rejected list which they have identified as confused, galaxy, 
PNe or foreground/background star. Noticeably 7 of their rejected sources do not have a match in 
new WISE within 1 ", and 4 do not have a match within 2" , therefore we cannot comment about these 
sources. 


3.3. Comparing with color-cuts 


We compared our candidates with the color-cut spaces described by Koenig et al. (2012). The 


color cuts are defined in multiple color spaces based on known different types of YSOs and other 
contaminant objects such as AGNs, galaxies with PAH-feature emissions and shock emission blobs. 


Koenig et al. ( 

2012 

) selections are employing the results from Spider surveys to identify various objects 

with color excess (Gutermuth et al. 

2009 Rebull et al. 2010 

) based on where they statistically lie on 


WISE + 2MASS color-color spaces. Similar to a selection technique, color-cut criteria in general can 
be adjusted to extract weaker disks with smaller IR excess or exclude them. In this section we compare 


our results with criteria defined by Koenig et al. (2012) which is frequently adopted in similar studies. 


Following Koenig et al. (2012), no AGNs were identified in our sample of candidates, but three 
sources (J03504369+3507088, J04022740+3057450, J03504333+3346024) might be considered as PAH 
emitting/star forming galaxies. These objects have no match in the literature. They are out of the 
Spitzer field and not covered by SDSS, therefore their identity cannot be determined by the spectrum. 
It worths to mention that Koenig color-cuts are defined based on the WISE preliminary data, while 
photometry of 25% of sources have had slight changes in all sky release. That may partially affect 
their selection criteria at the borders. 


Koenig et al. (2012) also used color-cuts to categorize their YSOs into different class types. 


Muench et al. (2007) (appendix A) explored the influence of dust extinction on a detection technique 
and presented that IRAC SED slope for a diskless K0 star requires Ay > 100 to inflect into a positive 
slope and even Ay > 200 cannot inflect <^ 5 . 8-24 of background stars into positive. Such large column 
densities within typical molecular clouds occur only in protosetllar envelop of embedded YSOs. But 
smaller values of dust extinction also may affect a at the classification of ClassI/ClassII candidates. In 
contrast mid-IR color cut techniques can be defined such that they minimize extinction-induced bias. 
The goal of this study is only identifying YSOs and not classifying them. We present the classification 
only to compare the slope-fitting and the color cut methods and will not report the class types for our 
final candidates. 


According to Muench et al. (2007), o^iraC3-8 — —2.66 corresponds to the predicted slope of 


M0 star photosphere. Accepting a typical precision of SED fits, they characterized sources with 
—2.56 < a < —1.80 as “anemic disks”. Transitions disks, disks with inner holes and heavily depleted 

































optically thin disks may locate in this category. Objects with aiRAC3-8 < —2.56 are considered stars, 
and objects with aiRAC3- 8 > —1-80 are accepted as IR excess objects (—1.80 < &irac3-8 < —0.5 as 
Class II and —0.5 < &irac3-8 as Class I). 


Figure [T] compares the disk classification by Muench et al. (2007) based on the slopes (a) and 
Koenig color cuts in one plot. Grey dots show all 48,692 WISE sources in our field with SNR>7 
for first three WISE bands. Black dots are a— stars. Yellow dots show a —anemic disks. Blue dots 
present Class II and red dots present Class I <a-disks. Based on Koenig classifications, green open 
circles present Class II and black open squares present Class I sources (selected in [3.4] — [4.6] and 
[4.6] — [12] color-color diagrams). There is a good agreement in disk identification type between two 
methods. However, color-cuts select less Class I protostars and many a —Class I sources are identified 
as Class II by color-cuts. 

We should remind that a selection fails to detect faint sources (Jmag >12) even if they have 
noticeable color excess. On the other hand color-cut method misses the weak or evolved disks, even if 
they are bright sources. For example anemic disks have rarely been distinguished by color-cuts. The 
comparison between the selected YSOs in Rebull et al. (2011) and this work in Figure [6] also presents 
how the faint sources even with strong disks have been missed in a selection, and how the bright 
sources with weak disks have been missed by color-cut selection. We suggest that to have a complete 
sample of faint to bright and weak to strong disks, both methods needed to be considered. 

To summarize, we establish our working sample of candidate young stellar objects as selected by 
the spectral energy distribution power-law slope in range of 3.4—12 /im from WISE all sky point source 
catalog. Using this method we can identify the weak and more evolved circumstellar disks if they are 
hosted by a bright protostar. These type of disks are barely identified with color-color diagrams. In 
contrast a^-u selection miss fainter sources, in particular those with J mag >12 sometimes even if they 
have strong disks. Color-color diagrams are capable to select many fainter sources but the results are 
contaminated (^ 67% for example, in Rebull’s primary selected sources) with background galaxies 
and other fake detections compared to a selection method. For example in pool of 1014 potential 


YSOs identified by Rebull et al. (2011), 686 were rejected as known galaxies or objects to be likely 


background galaxies in their final list. 


4. CANDIDATE ANALYSIS 


4.1. Known YSOs 


We cross matched our candidates with the SIMBAD database to identify known YSOs and pos¬ 
sible contaminants. Among 354 candidates passed our manual inspection as resolved point sources, 
J035523.11+310245.0 is a known massive X-ray binary and was removed. 12 other sources are known 
AGBs, Mira variables or other type of known evolved stars. 76 have matches with previously known 
YSOs or T-tau stars while 80 others are listed as YSO or T-tau candidates in SIMBAD. All these 
sources are categorized as known YSOs in our final catalog and listed in Table [3j 185 sources remain 
in our list as new YSO candidates, but 66 of them identified to be more likely evolved dusty stars 
and were separated from final new YSO candidates list (Table |4j We will discuss them in section 4.2 
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Finally 119 sources survived as new YSO candidates in our final catalog and they are listed in Table 
[2} The errors (range of Lower limit and upper limit for each parameters) are provided in online tables. 


4.2. Dusty evolved stars 


We checked the distribution of our sample in various color and magnitude diagrams to look 
for possible category and grouping objects. In particular 22 fi m is the best indicator to probe the 
embedded objects. Figure [8] shows 2MASS K mag plotted vs. WISE 22 /m 1 for all 353 point sources 
(J035523.11+310245.0, the known massive X-ray binary is removed). It was noticeable that sources 
are divided into two main populations. Sources with the same W4 magnitude are divided into two 
brightness branches in ¥^ mag . The lower group are usually bright stars with smaller slopes but they have 
statistically significant excess emission in longer wavelengths compared to main sequence stars. Known 
YSOs and YSO candidates lie in the upper group, with larger K magnitudes. Extreme brightness in 
K mag of the lower population suggests that they can be dusty evolved stars. The majority of them 
are out of Perseus cloud and have low extinctions. 


To characterize this “lower population” we compared them with the known AGB stars in our 
sample and in the literature. Grey dots in Figure [8] present all WISE sources in our survey field with 
SNR>7 for all four bands which mostly contain main sequence stars. Then we matched all known 
evolved stars in the Galaxy (including AGB stars, carbon stars and Mira variables) from SIMBAD 
with WISE catalog. Cyan dots show evolved stars which have a WISE match with SNR>10 in all 
four bands and we call it AGB branch in our diagram. The AGB branch in this plot is above the main 
sequence and perfectly matches with the “lower population”. 10 of the known evolved stars in our 
sample also lie in this area. We selected 66 objects from the “lower population” in K—W4 diagram as 
dusty evolved star candidates. These objects are shown by black crosses in Figure [8] left panel. Right 
panel in Figure [8] presents the K and M type stars in our sample which had known spectral type in the 
literature. While the majority of them locate within known and candidate YSOs, a smaller population 
lie within the AGB branch, indicating these objects are most likely evolved dusty stars rather than 
YSOs. 


In addition we checked our AGB candidate with Blum et al. (2006) color-magnitude diagrams 


from their Large Magellanic Cloud Spitzer survey. In Figures 3 — 6 of their paper they present the 
location of carbon stars and extreme AGB stars on [3.6] vs. J—[3.6], [8] vs. J—[8], [3.6] vs. [3.6] —[8] 
and [24] vs. [8] —[24]. We replaced [3.6] by W1 (3.4 fi m), [8] by W3 (12 fi m), and [24] by W4 (22 (i m) 
and applied their criteria. Although the wave bands are slightly different, our 66 AGB candidates 
reasonably matches the AGB stars criteria on their various plots. These candidates are listed in Table 

SI 

We also did a statistical estimation of how many AGB stars we expect to identify in our field. 


A recent WISE study (Tu & Wang 2013) estimated 470,000 of AGB stars in the Galaxy with con¬ 


tamination uncertainty of 20% of other sources including YSOs. Jackson et al. (2002) estimated a 
total of 200,000 AGB stars in the Galaxy. Using their volume density around solar neighborhood and 
integrating over 10 kpc (distance at which we can detect bright AGBs in our sample) we expect to 
have 100-400 AGB stars in our field. Therefore 66 likely AGB candidates plus 12 known evolved stars 
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in our field seems a reasonable number. 


4.3. Spatial Distribution 

Figure [9] presents the location of known and new YSO candidates, known and new candidate AGB 
stars and OB stars in 12° x 12° Per OB2 field. All the known YSOs (blue triangles) and the majority 
of new YSO candidates (red circles) are located within the nebulous region, bright in green color in 
Figure [5] New YSO candidates in North-East corner of the field follow the shell like dust feature and 
a large number of new YSOs in North-West side, toward the California nebula, are located where the 
remaining of the original Perseus cloud forms dense filamentary structures. 

OB stars and known AGBs are shown in yellow and black dots respectively. Black crosses present 
our AGB candidates discussed above. These candidates are expected to be foreground/background 
sources, therefore we expect them to be randomly distributed in the field. Figure [9] does not confirm 
a random distribution. AGB candidates shown with crosses are following the cloud structure and 
also concentrated in North-West. Some of AGBs also might be below our locus and not selected in 
our sample. Therefore we suggest that a fraction of these objects might be YSOs with inaccurate 
photometry or at the border of contamination with AGBs. 


5. MODELING SPECTRAL ENERGY DISTRIBUTION 

A common procedure to find the spectral energy distribution model that better reproduces data 
is using x 2 minimization: comparing a grid of models to the data points, finding the y 2 value using 
the measurement errors, and picking the model that minimizes the y 2 (e.g. Robitaille et al. (2007)). 
We call this solution the “Best Fit” model. However, besides the Best Fit solution, there might be 
many other models that have slightly larger y 2 but that are also a reasonable solution for our fitting 
problem. If we consider the distribution of various physical parameters from all models with y 2 smaller 
than a certain threshold, the peak of such distribution may not occur at the value given by the Best 
Fit. In other words, the Best Fit solution might not be located in the region of the parameter space 
where more models are a good representation of our data. Alternatively, we can assume that the 
model parameters are random variables, then derive probability distribution functions (PDFs) for 
them given our data, and find the parameter values that maximize those PDFs. In this paper we call 
the latter solutions the “Peak” values for the parameters, and we will use a Bayesian method to find 
them. 


Robitaille et al. 


(2007 


5.1. Bayesian approach 

In Bayesian inference, the model parameters are considered as random variables for which poste¬ 
rior probability distributions can be derived from the likelihood of the solutions (that we obtain from 
the data and the measurement errors) and the parameter priors, where we encode the prior knowledge 
in these parameters before any data has been taken. Once we derive the posterior PDFs, we can draw 
samples from them in order to visualize the solution and find the absolute PDF maximum for each 
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parameter. An efficient way to perform this sample is by using a Markov Chain Monte Carlo (MCMC) 
method that randomly steps across the parameter space and at each iteration decides whether to ac¬ 
cept the step based on an acceptance probability that depends on the ratio of probabilities between the 
current position and the proposed new position. This approach is particularly useful for multivariate 
problems like the one at hand here. 


Our Bayesian algorithm is very similar to the implementation of the Chiburst code for fitting 
SEDs of star-forming systems (Martfnez-Galarza et al. 2014). Here we will describe only the aspects of 


the code that are relevant for the present work. As a model grid, we use the SED models of Robitaille 


et al. (2006), that have calculated the radiative transfer for the emission from the YSO as it traverses 


the disk and/or envelope that surrounds it. Synthetic SEDs are produced for a broad range of stellar 
masses (0.1 — 5OM 0 ), ages (10 3 -10 7 yr), as well as disk and envelope sizes and geometries. Given 
a set of photometry and assuming this model grid is a fair description of the actual distribution of 
physical parameters in YSOs, we attempt to find the probability that a particular model represents 
the properties of the observed YSO. 


This allows us to explore the inherent degeneracies that arise when we attempt to fit a limited 
amount of photometric data with a multivariate model. By solving for the PDF rather than just 
finding solutions listed by increasing y 2 , we are able to visualize all the likely solutions for a particular 
object photometry, given the observational uncertainties and assuming that the parameter space of 
the model grid is a fair universal sample of the population of YSOs in the Perseus cloud. We focus 
on the determination of three main model parameters, namely the YSO stellar mass (m*), age (£*), 
and Ay in the line of sight towards our objects. To some extent, all other model parameters are 
determined by the selection of mass and age, or difficult to constrain without unavailable photometric 
bands. 


5.1.1. The Probability Distribution Functions 

Suppose that we have obtained photometry data (D) of a YSO at different bands, with certain 
observational uncertainties associated. Given our data, we can calculate P(0|D), the posterior PDF 
for the set of parameters 0 as the product of the likelihood P(D |0) and the prior P(M). In the case of 
normally distributed measurement errors, the likelihood P(D |0) can be obtained from the distribution 
of reduced y 2 values: 


P(D |0) = ^exp(—l/2 Xr 2 ed ) (1) 

where the sum is marginalized for each model parameter over all possible models with a given 
value of the parameter. The prior P(M) is a measure of any previous knowledge that we have on a 
particular parameter or set of parameters. For example, if we have reliable extinction measurements in 
the line of sight towards a particular YSO, then we can constrain the possible solutions to our problem 
by constructing a prior on ( Ay ) that is compatible with those extinction measurements. Finally, we 
need to apply a normalization factor to our posterior PDF to guarantee that the probability of at 
least one model being a representation for our YSO equals one. 
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5.1.2. Priors 


We are mostly interested in obtaining a distribution of masses, evolutionary stages and optical 
extinctions for a sample of newly identified YSOs in the Perseus region. For the mass and age of 
the YSOs we have adopted the flat, uniform priors already set by the model grid, with boundaries 
set by O.1M 0 < M* < 50 M 0 and 10 3 yr < t* < 10 7 yr. We also chose a flat, uniform prior for the 
disk inclination, since we expect most of the inclination effects to cancel out for a sample of many 
YSO models. We must emphasize here that these priors are discreet (i.e., not all possible values of 
mass, age and inclination are possible), but we consider that the sampling is such that any uncertainty 
due to sampling is smaller than the uncertainties imposed by the observational errors and the model 
degeneracies. Most of the other priors (disk and envelope mass, accretion rate, etc.) are set by the 
mass and age of a particular YSO, and have been thoroughly explained in (Robitaille et al. 2006). We 
adopt a distance of 320 kpc to the Perseus molecular cloud, but allow for a random scaling of the model 
YSOs to account for a distance uncertainty of 25%. We also obscure the models and characterize this 
obscuration by an Ay value, to account for any foreground extinction in the line of sight towards the 


Perseus cloud. We construct the prior on Ay using an extinction map produced by Lombardi et al. 


(2010). The prior is simple: we restrict the possible values of Ay to be within 0.4 dex of the value 


obtained from the extinction map. We will latter modify this prior to show the effect of the Ay on 
determining the mass of the YSOs. 


5.1.3. Stepping across the parameter space 

Calculating the posterior PDF at every single point of the model grid and for all possible values 
of distance scaling and Ay would be computationally time consuming, specially as more data points 
are added and additional model parameters (i.e. multiplicity of YSOs) are considered. Instead, we use 
a Monte Carlo Markov Chain (MCMC) approach to draw samples from the posterior PDF across the 
multidimensional space of parameters 0. We use the Metropolis algorithm with a uniform transition 
probability for the Markov Chain (i.e., given our current position in the parameter space, it is equally 
probable moving in any direction). After enough iterations of the MCMC, the histogram of model 
parameters in our chain should be a good representation of the posterior PDF. It is from this final 
histogram that we select the “Peak” solutions for our parameters. 


5.2. Examples 

Figure [TO] shows the best fits and associated posterior PDFs for several types of SEDs found in 
Perseus sample. The method is quite successful at finding solutions for almost all types of SEDs. More 
importantly, the plotted PDFs show the usefulness of the method in determining the uniqueness of 
the solution. 

It is important to realize the meaning of the best fit value as compared to the most likely solution 
as represented by the peak of the posterior PDF. They do not always coincide, which should not be 
surprising. The reason is that the best fit can be located in a region of the parameter space where 
not too many models agree with the data. The value obtained from the peak of the PDF, on the 
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other hand, provides the most likely value for a given parameter given all possible combinations of 
a parameter value with all the other parameters, constrained by the information contained in the 
priors. In a sense, the posterior PDF is a full description of all possible solutions. This is a significant 
improvement over the listing of the ten best fitting models, since these can be all located very near 
the best fit value, but not necessarily include the peak of the PDF. Furthermore, as we will see later, 
bimodal solutions are also possible due to degeneracies between the parameters, and usually the 10 
best fitting models all fall near one of the two possible solutions to the problem. 


5.2.1. Degeneracies 


Degeneracy is a common issue of all data-fitting problems that more than one combination of 
model parameters give a reasonable solution within the observational uncertainties. This problem of 
model degeneracies has been particularly neglected in the case of YSO SED fitting, partly because of 
the lack of enough data to break those degeneracies, but also due to a blind faith on the “best fit” 
solution, even for objects with only a few photometric observations. Our method allows to clearly 
visualize the degeneracies between model parameters and hence it helps us consider other possible 
solutions that might be hidden in the complexity of our multi-parametrical model, and that might be 
in better agreement with independent determinations of the parameters. It also allows us to explore 
how those degeneracies behave as a function of the priors, i.e., how additional information or additional 
data points can break these degeneracies. 


Figures [TT] and 12 shows an example of model degeneracies. Shown are two different realizations 
of the fit to the source J03284618+3116385. Two groups of solutions are possible, one with low 
foreground extinction (Ay < 1) and a low stellar mass, and one with higher Ay and higher mass. 
Physically, this can be understood as a need for a higher stellar mass to produce the same UV and 
optical flux when more obscuration is at work. In the unconstrained case shown in the figure, both 
fits are equally satisfying, but the solution with small Ay appears as more likely. Notice that the PDF 
for age remains unchanged. This implies that there is a strong degeneracy between stellar mass and 
extinction when fitting multi-wavelength SEDs of YSOs. Without a better prior for the Ay (or the 
stellar mass), we might choose the low mass-solution, but additional evidence could make the other 
solution more likely. For example, if we know that the Ay has to be greater than 1, then the low mass 
solution has to be neglected. This degeneracy is observed in many of the SEDs, in particular those of 
class Ha YSOs, as classified by (Gutermuth et al. 2009). Figure 13 shows the same effect much more 
clearly in the 2D m*-Ay plane. The choice of a particular Ay determines the estimated mass of the 
YSO. 


In fact, rather than thinking in terms of a “best fit”, we should think in terms of likelihood of a 
solution given certain constrains on the values and all possible combinations of parameters. 

It is possible to modify the prior on Ay by considering independent determinations of the physical 


conditions. We have obtained extinction values from the maps produced by Lombardi et al. (2010) 
and constrained the Ay prior accordingly. Since the extinction map have their own uncertainties, we 
allow our Ay values to vary within 0.4 dex of the value obtained from the maps. By doing so, we 
hope to break the degeneracy observed in Figure [13| for at least some of the objects. It is important 
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to notice here that the values of Ay measured from the map represent in fact only an upper limit to 
the total extinction, since part of the measured Ay can actually be coming from behind the source 
of interest. By selecting a range of values around the measured Ay , we are assuming that most of 
the extinction we see is due to foreground material between us and the source. Figure [T3] shows the 
solution and the resulting marginalized posterior PDFs after the prior has been modified as described, 
for the same object of Figures [TT] and [l2| As expected, the additional information of the prior breaks 
the existing degeneracy and leaves only the most massive of the two solutions. This solution is in fact 
in better agreement with the expected mass for YSOs that we would be able to detect at the distance 
to the Perseus cloud. Of course, another way to break degeneracies might be by including additional 
data-points at longer wavelengths that are compatible with only one of the two possible solutions. 
Photometry from the PACS instrument onboard the Herschel Space Observatory will be instrumental 
for this task. 


5.3. Results and discussion 


5.3.1. SED slope analysis 


We have applied our SED fitting method to 275 identified YSOs (total of known and new candi¬ 
dates) in order to perform a census of masses and evolutionary stages. The calculated parameters are 
highly dependent on the extinction, and therefore we had to estimate the dust absorptions in order 
to constrain Ay before the SED fitting was applied. The extinction values are calculated based on 
2MASS K band extinction map provided by Lombardi et al. (2010). Figure 15 presents the location 
of high (Ay > 2) and low (Ay < 2) extinction regions. Currently active star forming regions such as 
NGC 1333, IC 348, L1448, L1455 and Taurus all are embedded within Perseus molecular clouds and 
have higher extinctions. 


Figures [16] presents the distribution of calculated 0 ^ 3 . 4-12 and position of candidates in three 
distinct groups. For the first group, 0 ^ 3 . 4-12 peaks at ~ —2 and these sources are noted by blue circles 
in right panel for objects with 0 ^ 3 . 4-12 < —1.75. This group have more evolved disks and are spread 
over the field with a concentration at North-West of the field, above the California nebula. 


The main peak at 0 ^ 3 . 4-12 ~ —1 contains many YSOs with stronger disks. Green triangles present 
this group, selected as —1.75 < <^ 3 . 4-12 < 0 in the right panel. They are more concentrated within the 
Perseus and Taurus star forming regions and molecular cloud. 


Candidates with very strong disks are presented with red dots and are selected as 0 < 0 ^ 3 . 4 - 12 - 
They are highly concentrated within active star forming clusters, NGC1333, IC348, and Taurus with 
a few spread over the field. It is not clear if they are formed individually or escaped from their original 
birthplace. 


Distinct population in Perseus region has been noted and discussed in several studies. For example 


Herbig (1998) found evidence of different generation of stars within and around IC348 cluster. Also 


suggested by a hipparcos study (de Zeeuw et al. (1999) and references therein), and later confirmed by 


Belikov et al. (2002) Per OB2 contains two kinematically separated subgroups. Belikov et al. (2002) 


identified more than 800 members for Per OB2 within a 50 pc diameter. We also have identified three 
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separated evolutionary stages which will be discussed more in next section. 


5.3.2. Parameters from SED fitting 

Using the Bayesian MCMC method described, we have obtained best fit and peak values for all 
the YSOs in our sample. We focus on the peak values for stellar mass, age, total luminosity and 
extinction. Given the limited amount of information at wavelengths longer than 22 /im, we do not 
attempt to constrain the disk/envelope properties with accurate details. Instead, we are interested 
in studying the distribution of masses and ages in our wide select field that includes NGC 1333 and 
IC 348, as well as another less dense concentration of YSOs located near California nebula. 


We should note that the derived parameters using Robitaille et al. (2007) models are “model- 


dependent-fit parameters” and not physical values. The YSO grid of models are randomly sampled 


in mass and age, and then use evolutionary tracks from Siess et al. (2000) to calculate all the other 
stellar properties. This means that the “age” derived here is not measured with respect to any physical 
stellar process, such as the start of deuterium burning. Rather, the reported ages are indications of the 
relative evolutionary stages of the members which are dependent on model assumptions and priors. 


Willis et al. (2013) found a noticeable difference in total mass of their large complete sample of 


YSOs in NGC6334, derived from the same set of models used in this paper, compared to the expected 
total mass for a complete Kroupa IMF. Robitaille models employed in this study are simulated based 
on observation and physical parameters, therefor results are valid as a relative ratio of the evolutionary 
stage of all sources in the survey that satisfies the goal of this study. 


Figure [14| shows the age probability vs. mass probability for 2MASS J03311069+3049405. The 
blue cross which presents the best fit is very close to the most probable models for both mass and 
age in left panel. In some cases, however, the best fit model is very different from the probability 
peak. 2MASS 03253790+3108207 (right panel) is a good example where Best Fit and Peak Value 
do not match well. There are two main reasons for this mismatch: first, the PDFs plotted here are 
marginalized versions of the PDFs, and are therefore integrated over all other parameters not shown 
in the plot. It is therefore likely that the absolute maximum of the multi-dimensional PDF does 
not always coincide with the maximum of a marginal distribution. Also, degeneracies in the models, 
specially given the lack of far-IR data points, might produce multiple maxima in the probability 
distributions. 

We also realized that the correct estimation of the optical extinction is very important to constrain 
the mass and other physical parameters. The fitting tool lets Ay to accept the values in a wide range, 
to find the best SED models. But not constraining Ay will dramatically affect the final results. We 
fixed this problem by constraining Ay from observed values of extinctions toward each individual 


candidate using the Ak value from the extinction map provided by Lombardi et al. (2010) toward 


Taurus and Perseus star forming regions, based on 2MASS counts. We have recorded the best fit 
values for the parameters for each source, as well as the most likely solution or the peak value as 
obtained from the peak of the posterior PDF. 


Figure mi presents the distribution of estimated ages for our candidates and the age spatial 















-16- 


distribution in the field. Three distinct age groups are noticeable in left histogram: age<l Myr, 
l<age<5 Myr and 5 Myr<age. As expected, the younger YSOs are having stronger disks and lo¬ 
cate within the high extinction inner parts of the Perseus cloud and NGC1333 and IC348 clusters. 
Noticeably the North-Western association do not contain any young YSOs and only a few mid-age 
ones. 

The estimated masses for our YSO candidates varies in a wide range from 0.1 to 5 M 0 . Figure [l8| 
presents the mass distribution and the location of low mass and high mass candidates. The majority 
of YSOs with masses larger than 1 M 0 are located within areas with larger Ay, i.e. within Perseus 
molecular cloud and particularly star-forming clusters. The massive YSOs from older populations 
spread over the field are probably too evolved to have strong, detectable disks. 


6. summary 


We have performed a census of YSO candidates in Perseus OB2 association covering ^144 deg 2 
using WISE catalog. We derived physical characteristics of all identified YSOs within the region by 
employing other optical and infrared data including 2MASS, Spitzer , SDSS, PPMXL and APASS. 
Following Lada et al. (2006) Muench et al. (2007) we calculated the SED slope in range of 3 — 12 /im 


for 48,692 sources which had a S/N<7 in first three WISE bands. 669 point sources survived to 
have slopes larger than the 5 a above the peak in bins of 0.5 2MASS J magnitude. Removing known 
non-stellar extended objects such as galaxies and PNe we still had some sources remained in our 
sample that could not be identified as point sources in optical or infrared images or did not show 
point source profiles in at least three out of four WISE bands. Among 354 remained point sources 
one was identified as a massive X-ray binary and 12 as evolved stars such as carbon stars or Mira 
variables. 156 of remained sources were previously identified as YSO or YSO candidates. In this work 
we present 119 new YSO candidate toward the Perseus region. We also identified 66 new point sources 
with infrared excess emission but brighter than normal YSOs in the region. We separated them as 
likely AGB and evolved star candidates. 

The majority of known candidates are concentrated toward active star forming regions such as 
Taurus, IC348 and NGC1333 clusters. We add more YSO candidates in these regions which had poor 
photometry in c2d or other previous studies to be confirmed as YSOs. New candidates also follow the 
remaining gas of the original Perseus cloud with a concentration toward North-West near California 
nebula. In total, new candidates have weaker disks and are more scattered within the field, while 
previously known sources locate within high Ay regions in Perseus molecular cloud. 


We employed the SED fitting models described by Robitaille et al. (2007), but used a Markov 


Chain Monte Carlo method to explore the large parameter space of the model grid. Instead of selecting 
the model with the lowest y 2 solution as the best fit, we accept the physical parameters from the most 
probable fitted models. Similar to the slope (A 3 . 4 - 12 ) distribution, derived mass and ages also present 
separated populations. <^ 3 . 4-12 histogram shows two separated peaks at ~ —2 and ~ — 1 . Candidates 
with slopes larger than zero (i.e. very strong disks) are mainly found within IC348 and NGC1333 
clusters with a few scattered in the field. The mid-slope population which peaks around —1 are also 
located within the Perseus cloud and nebulous region. 
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As expected, younger sources (age<l Myr) are located within active star forming clusters while 
older population (age>5 Myr) are scattered within the field with a noticeable concentration toward 
North-West. In contrast both low mass (M<1 M 0 ) and higher mass (M>1 M 0 ) candidates are found 
equally scattered within the filed or located within the Perseus cloud. 

Finally, we compared our method of selecting circumstellar disks with other YSO selecting meth¬ 
ods based on optical and infrared color and magnitudes. While our a method is very successful to 
identify bright sources with weak disks that will be filtered in other methods, our method misses YSOs 
with strong disks if they are faint and in particular fainter than 3 ma g — 12. A combination of both a 
and color-cut selection would be a complete method to identify both strong and weak disks. 
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Fig. 1.— Surveyed Region and its coverage by a sample of surveys: WISE (blue), Spitzer (c2d) 
(yellow), SDSS (green), UKIDSS (red), Gallex (grey), XMM (black) and Chandra (magenta). 
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Fig. 2.— Field coverage by SDSS, APASS and WISE 
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Fig. 3.— IRAC « 3-8 vs. WISE 0 ^ 3 - 12 • Solid blue line presents a linear least square fit with slope of 
0.92 ± 0.2 (excluding three data points at right bottom corner with poor IRAC photometry) and red 
dashed line presents equal values. 
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2MASS Jmag 


Fig. 4.— a' 3 . 4 _i 2 versus 2MASS J magnitude were used to identify the YSO candidates. Grey dots 
show all WISE point sources with SNR>7 in first three bands. Black locus separates sources with 
Q' 3 . 4 _i 2 which are 5 a above the Gaussian peak of slope distribution in each J ma ^ 0.5 magnitude bins. 
New identified YSO candidates are presented with open red circles while known candidates from 
literature are shown in blue triangles. Known AGB stars are presented in black bold dots in this plot 
while AGB candidates identified in this study are presented in black crosses. Yellow dots present OB 
stars in Per OB2. 
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Fig. 5.— YSO candidates overlaid on multi band WISE image of Perseus region (4.6/im in blue, 
12/un in green and 22/rni in red). The brightness of each source presents its 0 ^ 3 . 4-12 slope. Larger 
slope presumably means stronger disks. Picture created using Microsoft Worldwide Telescope (http : 
/ / worldwidetelescope.org) 
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Fig. 6.— Comparing our results with Rebull et al. (2011) study of Taurus-Auriga star forming region. 
Black triangles present our sample of 354 candidates in Perseus. Green, blue and pink dots presents 
known YSO candidates, new YSO candidates and rejected candidates respectively from Rebull et.al. 
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Fig. 7.— Comparing disk classification based on the slopes (a) and Koenig color cuts in one plot. 
Black dots show a —stars (a < —2.56) and yellow dots show a —anemic disks (—2.56 < a < —1.8). 
Blue dots show Class II (—1.8 < a < —0.5) and red dots present Class I cr-disks (a > —0.5). Selected 
by [3.4] — [4.6] and [4.6] — [12] color-color diagrams, green open circles present Class II and black open 
squares present Class I sources. 
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Fig. 8.— Left: our candidate with infrared excess are divided into two groups. The brighter sources 
with brighter K mag match with known Galactic AGB stars (cyan dots) in this plot and are likely AGB 
stars and not YSOs. They are noted with black crosses. Right: M and K type stars in our sample 
similarly follow the same pattern, empowering the suggestion that sources noted by cross in left panel 
are likely M and K type evolved dusty stars. 















-26- 



Fig. 9.— Distribution of YSO, AGB and OB stars in Per OB2. As expiated YSOs are concentrated 
toward star forming clusters (IC348, NGC1333) and more numerous within the Perseus molecular 
cloud. AGB candidates are expected to be randomly distributed in the field but they are following 
the cloud structure and also concentrated in North and North-West of the field. 
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Fig. 10.— Multi-wavelength SEDs for a selection of objects in Perseus. The black diamonds are the 
extracted photometry, while the red line is the best fit (the one that minimizes the y 2 ). Also shown 
are the posterior PDFs for m*, £*, L tot and Ay, with the best fit values represented by the dotted 
lines. 
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Fig. 11.— The m*-Ay degeneracy. Two possible solutions for a single set of photometry. The PDFs 
show the two possible solutions, while the best fit in each case tends to select only one. 
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Fig. 12.— The m*-Ay plane with contours of normalized probability showing the degeneracy between 
these two model parameters. The blue cross indicates the best fit for this particular realization. 



Wavelength [/zm] 


Fig. 13.— Only one solution is possible when the Ay prior is modified to account for additional 
evidence on the extinction. 
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Fig. 14.— Mass-Age probability distribution for all fitted SED models for two different YSO candi¬ 
dates. The brightest region presents the most probable values from different fitted SED models. In 
the right panel the best-fit values (blue cross)is close to the most probable values but in the right 
panel there is a large difference. 



Fig. 15.— Spatial distribution of extinction for 275 YSO candidates. Ays are calculated based on 
2MASS K band extinction map by Lombardi et al. 
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Fig. 16.— Distribution of « 3 . 4 _i 2 for 275 YSO candidates. The location of members of three noticeable 
groups are shown in right panel, a < —1.75 in blue circles, — 1.75<<a<0in green triangles and 
0 < a in red dots. 




Age (Myr) 


RA 


Fig. 17.— Age distribution for 275 YSO candidates. The location of members of three noticeable 
groups are shown in right panel, age < 1 Myr in red dots, 1 < age < 5 Myr in green triangles and 
5 < age in blue circles. 
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Fig. 18.— Mass distribution for 275 YSO candidates. The location of two noticeable groups are 
shown in right panel. M< 1 M 0 in red dots and M> 1 M 0 in blue squares. 
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J033805.39+325428.5 54.5224599 32.9079195 -1.16 ± 0.07 0.62 0.47 ± 0.07 6.58 2.91 ± 0.3 

J033900.55+294145.7 54.7523308 29.6960465 V* V1185 Tau -0.08 ± 0.04 3.3 2.17 ± 0.1 4.35 7.95 ± 0.13 

J033909.97+322421.8 54.7915742 32.4060745 -2.15 ± 0.08 0.76 0.85 ± 0.26 3.51 4.86 ± 0.23 

J033941.31+361604.3 54.9221469 36.2678697 GSC 02367-01706 -1.0 ± 0.15 2.76 2.8 ± 0.0 0.1 0.11 ± 0.0 

J033953.88+320742.6 54.9745076 32.1285048 -0.53 ± 0.07 1.25 0.3 ± 0.16 4.39 8.39 ± 0.22 
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J035216.28+332422.1 58.0678615 33.4061622 HD 279128 -1.27 ± 0.05 2.59 2.89 ± 0.04 2.71 2.1 ± 0.28 

J035224.19+361221.4 58.1008192 36.2059517 HD 279075 -2.5 ± 0.05 1.91 1.77 ± 0.14 8.85 8.28 ± 0.25 

J035234.42+353358.9 58.1434471 35.5663787 -1.12 ± 0.04 0.36 0.35 ± 0.24 4.49 2.27 ± 0.2 

J035323.82+271838.3 58.3492848 27.310649 -0.99 ± 0.04 0.11 0.38 ± 0.14 7.59 6.91 ± 0.18 

J035357.09+350239.2 58.4878805 35.0442238 -1.13 ± 0.04 1.87 0.69 ± 0.03 3.61 4.17 ± 0.3 
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101 J040227.39+305745.3 60.6141306 30.9625879 0.39 ± 0.04 0.24 0.3 ± 0.13 9.42 9.22 ± 0.05 

102 J040259.96+315703.9 60.7498591 31.9510947 -2.05 ± 0.04 0.6 0.51 ± 0.11 0.31 2.23 ± 0.27 

103 J040323.26+360327.2 60.8469493 36.0575797 -2.33 ± 0.06 1.36 1.36 ± 0.08 9.49 4.97 ± 0.2 

104 J040401.78+271545.4 61.0074512 27.2626172 IRAS 04010+2707 -1.92 ± 0.05 0.81 0.78 ± 0.3 0.17 0.28 ± 0.21 

105 J040506.06+361129.2 61.2752807 36.1914563 -1.63 ± 0.07 0.95 0.89 ± 0.13 7.92 6.9 ± 0.13 
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J032921.87+311536.2 52.3411571 31.2600729 EM* LkHA 271 -1.53 ± 0.04 1.08 1.82 ± 0.09 7.46 2.46 ± 0.28 

J032923.15+312030.3 52.3464585 31.3417503 EM* LkHA 355 -0.89 ± 0.04 1.48 0.96 ± 0.32 8.59 4.59 ± 0.19 

J032923.23+312653.0 52.3468186 31.4480694 2MASS J03292322+3126531 -0.53 ± 0.04 0.39 0.39 ± 0.1 5.14 3.8 ± 0.38 

J032925.92+312640.1 52.3580264 31.4444756 2MASS J03292591+3126401 -0.78 ± 0.03 1.48 1.87 ± 0.09 2.24 2.38 ± 0.4 

J032926.79+312647.4 52.3616596 31.44652 2MASS J03292681+3126475 -1.94 ± 0.05 1.6 1.39 ± 0.16 5.2 3.44 ± 0.25 
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J034114.12+315946.1 55.308873 31.996148 2MASS J03411412+3159462 -1.88 ± 0.05 1.97 2.04 ± 0.14 9.79 8.31 ± 0.1 

J034157.45+314836.6 55.4894069 31.810175 [EDJ2009] 314 -1.27 ± 0.04 2.28 2.21 ± 0.36 7.82 8.39 ± 0.62 

J034157.76+314800.7 55.490705 31.8002143 [EDJ2009] 315 -0.89 ± 0.04 0.42 1.73 ± 0.1 0.45 0.87 ± 0.37 

J034204.35+314711.4 55.5181477 31.7865217 [EDJ2009] 319 -0.84 ± 0.04 0.43 0.21 ± 0.39 1.81 1.83 ± 0.44 

J034219.28+314326.8 55.5803535 31.7241352 2MASS J03421927+3143269 -0.8 ± 0.04 1.08 1.17 ± 0.26 7.46 0.98 ± 0.34 



Order Catalog Number RA Dec Name a Mass bf Mass pea k Ag qbf Age pea k 

(deg) (deg) (M 0 ) (M 0 ) (Myr) (Myr) 
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136 J034516.34+320619.9 56.3180861 32.1055416 EM* LkHA 99 -0.8 ± 0.04 2.29 1.5 ± 0.15 7.17 0.59 ± 0.17 

137 J034520.45+320634.4 56.3352286 32.1095804 2MASS J03452046+3206344 -1.15 ± 0.04 2.16 2.05 ± 0.24 2.02 2.44 ± 0.39 

138 J034525.14+320930.3 56.3547711 32.1584269 2MASS J03452514+3209301 -0.79 ± 0.04 0.57 0.53 ± 0.27 3.11 1.11 ± 0.3 

139 J034536.85+322556.8 56.4035733 32.4324705 EM* LkHA 329 -1.31 ± 0.04 1.25 1.43 ± 0.14 0.65 0.67 ± 0.28 

140 J034548.28+322411.8 56.4512063 32.4032995 IRAS 03426+3214 -1.82 ± 0.06 2.46 3.11 ± 0.1 3.92 0.5 ± 0.31 
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